Explorar os dados sobre os mananciais utilizados para abastecimento público na Região Metropolitana de São Paulo.
Perguntas:
Qual é o sistema mais importante para o abastecimento público atualmente?
Atualmente, qual é o sistema com volume armazenado mais baixo? e mais alto?
Considerando os gráficos gerados, você diria que o Cantareira se recuperou da crise hídrica de 2014-2015?
Fonte dos dados:
library(tidyverse)
library(knitr)
library(lubridate)
library(janitor)
library(broom)
library(gghighlight)
library(skimr)
library(DT)
# ler o arquivo
mananciais <-
read_delim(
"https://github.com/beatrizmilz/mananciais/raw/master/inst/extdata/mananciais.csv",
delim = ";",
escape_double = FALSE,
col_types = cols(data = col_date(format = "%Y-%m-%d")),
locale = locale(decimal_mark = ",", grouping_mark = "."),
trim_ws = TRUE
)
Quais colunas estão disponíveis?
glimpse(mananciais)
## Rows: 50,997
## Columns: 8
## $ data <date> 2022-07-12, 2022-07-12, 2022-07-12, 2022-07-12, 2…
## $ sistema <chr> "Cantareira", "Alto Tietê", "Guarapiranga", "Cotia…
## $ volume_porcentagem <dbl> 38.4, 58.1, 72.1, 79.0, 95.5, 43.0, 85.3, 38.5, 58…
## $ volume_variacao <dbl> -0.1, -0.2, -0.4, -0.2, -0.1, -0.2, -0.3, -0.1, -0…
## $ volume_operacional <dbl> 376.97898, 325.41513, 123.36649, 13.03291, 107.103…
## $ pluviometria_dia <dbl> 0.0, 0.0, 0.0, 0.2, 0.2, 0.2, 0.2, 0.1, 0.2, 0.0, …
## $ pluviometria_mensal <dbl> 0.3, 1.9, 0.2, 1.2, 0.6, 2.8, 1.2, 0.3, 1.9, 0.2, …
## $ pluviometria_hist <dbl> 46.5, 46.9, 41.5, 51.3, 54.2, 91.3, 77.7, 46.5, 46…
skim(mananciais)
| Name | mananciais |
| Number of rows | 50997 |
| Number of columns | 8 |
| _______________________ | |
| Column type frequency: | |
| character | 1 |
| Date | 1 |
| numeric | 6 |
| ________________________ | |
| Group variables | None |
Variable type: character
| skim_variable | n_missing | complete_rate | min | max | empty | n_unique | whitespace |
|---|---|---|---|---|---|---|---|
| sistema | 0 | 1 | 5 | 12 | 0 | 7 | 0 |
Variable type: Date
| skim_variable | n_missing | complete_rate | min | max | median | n_unique |
|---|---|---|---|---|---|---|
| data | 0 | 1 | 2000-01-01 | 2022-07-12 | 2011-08-20 | 8229 |
Variable type: numeric
| skim_variable | n_missing | complete_rate | mean | sd | p0 | p25 | p50 | p75 | p100 | hist |
|---|---|---|---|---|---|---|---|---|---|---|
| volume_porcentagem | 0 | 1 | 66.57 | 25.72 | -24.30 | 47.70 | 69.20 | 88.7 | 110.50 | ▁▂▇▇▇ |
| volume_variacao | 0 | 1 | 0.00 | 0.59 | -6.40 | -0.20 | -0.10 | 0.1 | 16.60 | ▁▇▁▁▁ |
| volume_operacional | 0 | 1 | 151.27 | 185.47 | -238.84 | 13.96 | 97.02 | 200.8 | 987.05 | ▁▇▂▁▁ |
| pluviometria_dia | 0 | 1 | 4.05 | 9.86 | 0.00 | 0.00 | 0.20 | 2.8 | 222.50 | ▇▁▁▁▁ |
| pluviometria_mensal | 0 | 1 | 64.14 | 71.03 | 0.00 | 9.80 | 41.20 | 94.4 | 557.60 | ▇▂▁▁▁ |
| pluviometria_hist | 245 | 1 | 132.26 | 68.40 | 23.50 | 75.00 | 123.40 | 189.1 | 298.90 | ▇▇▅▅▂ |
Criar a coluna ano, mes_ano, e ordenar a coluna de sistema:
mananciais <- mananciais |>
mutate(
ano = year(data),
mes_ano = floor_date(data, "month"),
sistema = fct_relevel(
sistema,
levels = c(
"Cantareira",
"Alto Tietê",
"Guarapiranga",
"Rio Grande" ,
"São Lourenço" ,
"Cotia" ,
"Rio Claro"
)
)
)
Quais são os sistemas existentes na base? Quantas observações temos para cada Sistema?
mananciais |>
count(sistema) |>
kable()
| sistema | n |
|---|---|
| Cantareira | 8229 |
| Alto Tietê | 8229 |
| Guarapiranga | 8229 |
| Rio Grande | 8229 |
| São Lourenço | 1623 |
| Cotia | 8229 |
| Rio Claro | 8229 |
mananciais |>
group_by(sistema) |>
filter(data == min(data)) |>
select(sistema, data) |>
kable()
| sistema | data |
|---|---|
| São Lourenço | 2018-02-01 |
| Cantareira | 2000-01-01 |
| Alto Tietê | 2000-01-01 |
| Guarapiranga | 2000-01-01 |
| Cotia | 2000-01-01 |
| Rio Grande | 2000-01-01 |
| Rio Claro | 2000-01-01 |
Quais são os sistemas com o maior volume operacional armazenado?
mananciais |>
filter(data == max(data)) |>
arrange(desc(volume_operacional)) |>
select(sistema, volume_operacional) |>
kable()
| sistema | volume_operacional |
|---|---|
| Cantareira | 376.97898 |
| Alto Tietê | 325.41513 |
| Guarapiranga | 123.36649 |
| Rio Grande | 107.10315 |
| São Lourenço | 75.72723 |
| Cotia | 13.03291 |
| Rio Claro | 5.87214 |
mananciais |>
filter(sistema == "Cantareira") |>
ggplot() +
geom_line(aes(x = data, y = volume_porcentagem)) +
theme_minimal() +
labs(x = "Ano", y = "Volume (%)")
mananciais |>
filter(sistema == "Cantareira") |>
ggplot() +
geom_line(aes(x = data, y = volume_porcentagem)) +
theme_minimal() +
labs(x = "Ano", y = "Volume (%)") +
gghighlight(ano %in% c(2014:2015))
# gghighlight(mes_ano >= "2013-04-01", mes_ano < "2016-01-01")
Adapte o código abaixo para criar um gráfico para o sistema Alto Tietê:
mananciais |>
filter(sistema == "Cantareira") |>
ggplot() +
geom_line(aes(x = data, y = volume_porcentagem)) +
theme_minimal() +
labs(x = "Ano", y = "Volume (%)")
mananciais |>
summarise(min = min(volume_porcentagem),
max = max(volume_porcentagem),
media = round(mean(volume_porcentagem), 1),
variancia = round(var(volume_porcentagem), 1),
desvio_padrao = round(sd(volume_porcentagem), 1)) |>
kable()
| min | max | media | variancia | desvio_padrao |
|---|---|---|---|---|
| -24.3 | 110.5 | 66.6 | 661.8 | 25.7 |
mananciais |>
group_by(ano, sistema) |>
summarise(min = min(volume_porcentagem),
max = max(volume_porcentagem),
media = round(mean(volume_porcentagem), 1),
variancia = round(var(volume_porcentagem), 1),
desvio_padrao = round(sd(volume_porcentagem), 1)) |>
datatable()
Adapte o código a seguir para criar uma tabela de resumo para a variável de Volume Operacional:
mananciais |>
group_by(ano, sistema) |>
summarise(min = min(volume_porcentagem),
max = max(volume_porcentagem),
media = round(mean(volume_porcentagem), 1),
variancia = round(var(volume_porcentagem), 1),
desvio_padrao = round(sd(volume_porcentagem), 1)) |>
datatable()
modelo_linear <- lm(volume_variacao ~ pluviometria_dia, data = mananciais)
modelo_linear
##
## Call:
## lm(formula = volume_variacao ~ pluviometria_dia, data = mananciais)
##
## Coefficients:
## (Intercept) pluviometria_dia
## -0.13892 0.03504
tidy(modelo_linear) |> kable()
| term | estimate | std.error | statistic | p.value |
|---|---|---|---|---|
| (Intercept) | -0.1389243 | 0.0022709 | -61.17589 | 0 |
| pluviometria_dia | 0.0350441 | 0.0002131 | 164.41319 | 0 |
glance(modelo_linear) |> kable()
| r.squared | adj.r.squared | sigma | statistic | p.value | df | logLik | AIC | BIC | deviance | df.residual | nobs |
|---|---|---|---|---|---|---|---|---|---|---|---|
| 0.3464416 | 0.3464288 | 0.4743759 | 27031.7 | 0 | 1 | -34329.33 | 68664.66 | 68691.18 | 11475.53 | 50995 | 50997 |
mananciais |>
group_by(sistema) |>
summarise(
modelo = lm(volume_variacao ~ pluviometria_dia, data = cur_data()) |>
tidy() |>
select(term, estimate) |>
pivot_wider(names_from = term, values_from = estimate),
r2 = lm(volume_variacao ~ pluviometria_dia, data = cur_data()) |>
glance() |>
pull(adj.r.squared) |>
round(2)
) |>
unnest(cols = c(modelo)) |>
clean_names() |>
mutate(
modelo = str_glue(
"volume_variacao = {round(intercept,3)} + {round(pluviometria_dia,3)} * pluviometria_dia"
)
) |>
select(sistema, modelo, r2) |>
arrange(desc(r2)) |>
kable()
| sistema | modelo | r2 |
|---|---|---|
| Rio Claro | volume_variacao = -0.279 + 0.048 * pluviometria_dia | 0.51 |
| Cotia | volume_variacao = -0.125 + 0.04 * pluviometria_dia | 0.41 |
| Rio Grande | volume_variacao = -0.142 + 0.035 * pluviometria_dia | 0.38 |
| Guarapiranga | volume_variacao = -0.105 + 0.031 * pluviometria_dia | 0.29 |
| Alto Tietê | volume_variacao = -0.074 + 0.02 * pluviometria_dia | 0.28 |
| Cantareira | volume_variacao = -0.07 + 0.018 * pluviometria_dia | 0.21 |
| São Lourenço | volume_variacao = -0.103 + 0.035 * pluviometria_dia | 0.12 |
mananciais_r2 <- mananciais |>
group_by(sistema) |>
summarise(
r2 = lm(volume_variacao ~ pluviometria_dia, data = cur_data()) |>
glance() |>
pull(adj.r.squared) |>
round(2),
max_pluviometria_dia = max(pluviometria_dia),
max_volume_variacao = max(volume_variacao),
)
mananciais |>
ggplot(aes(y = volume_variacao, x = pluviometria_dia)) +
geom_point(aes()) +
geom_smooth(method = "lm") +
facet_wrap(vars(sistema)) +
geom_text(aes(
x = max(max_pluviometria_dia)*0.80,
y = max(max_volume_variacao)*0.9,
label = paste0("R2 = ", r2)
), data = mananciais_r2) +
theme_minimal() +
labs(x = "Pluviometria do dia", y = "Variação do volume diário")
A base original apresenta dados a partir do ano 2000. Altere o código abaixo para filtrar os anos pós crise hídrica (a partir de 2016).
Repare no R2 gerado para o Sistema Cantareira neste gráfico (pós crise) e no gráfico anterior: o que você acha que causa essa diferença no R2?
manancias_filtrado <- mananciais |>
filter(ano >= 2000)
mananciais_r2_filtrado <- manancias_filtrado |>
group_by(sistema) |>
summarise(
r2 = lm(volume_variacao ~ pluviometria_dia, data = cur_data()) |>
glance() |>
pull(adj.r.squared) |>
round(2),
max_pluviometria_dia = max(pluviometria_dia),
max_volume_variacao = max(volume_variacao),
)
manancias_filtrado |>
ggplot(aes(y = volume_variacao, x = pluviometria_dia)) +
geom_point(aes()) +
geom_smooth(method = "lm") +
facet_wrap(vars(sistema)) +
geom_text(aes(
x = max(max_pluviometria_dia)*0.80,
y = max(max_volume_variacao)*0.9,
label = paste0("R2 = ", r2)
), data = mananciais_r2_filtrado) +
theme_minimal() +
labs(x = "Pluviometria do dia", y = "Variação do volume diário")
https://abjur.github.io/livro/ - Livro em construção sobre Jurimetria